% 反向传播代码
function [W1, W2, W3, error, accuracy] = BP(W1, W2, W3, alpha, D, imageData, accuracy)
    %卷积
    img_conv1 = zeros(20, 20, 20);
    for k = 1:20
        img_conv1(:, :, k) = filter2(W1(:, :, k), imageData, 'valid'); 
    end
    %ReLU激活
    img_act = max(0, img_conv1);
    %池化
    img_pool = (img_act(1:2:end, 1:2:end, :) + img_act(2:2:end, 2:2:end, :) +     
    img_act(1:2:end, 2:2:end, :) +img_act(2:2:end, 1:2:end, :)) / 4;
    %将img_pool转换成列向量（2000*1)
    img_input = reshape(img_pool, [], 1);
    %第一个隐层的输出
    v1 = W2 * img_input;
    y1 = max(0, v1);
    %输出层的输出
    v2 = W3 * y1;
    y2 = softmax(v2);
    %观察模型是否训练准确
    [value1, index1] = max(y2);
    [value2, index2] = max(D);
    if index1 == index2
        accuracy = accuracy + 1;
    end
    % 计算交叉熵函数
    error = sum(- D .* log(y2) - (1 - D) .* log(1 - y2)) / 10;
    %误差反向传播过程
    %计算输出层的delta
    e2 = D - y2;
    delta2 = e2;
    %计算第一层隐藏层的delta1
    e1 = W3' * delta2;
    delta1= (y1 > 0) .* e1;
    %计算输入层(reshape层)的e
    e = W2' * delta1;
    %将输入层的误差进行reshape，以便于误差进一步反向传播穿过池化层和卷积层
    E2 = reshape(e, size(img_pool));
    %将池化层的误差传播到卷积层
    E1 = zeros(size(img_act));
    E2_4 = E2 / 4;
    E1(1:2:end, 1:2:end, :) = E2_4;
    E1(1:2:end, 2:2:end, :) = E2_4;
    E1(2:2:end, 1:2:end, :) = E2_4;
    E1(2:2:end, 2:2:end, :) = E2_4;
    delta = (img_act > 0) .* E1;
    dW1 = zeros(9, 9, 20);
    for k = 1:20
        dW1(:, :, k) = alpha * (filter2(delta(:, :, k), imageData, 'valid'));
    end
    %更改权重
    W1 = W1 + dW1;
    W2 = W2 + alpha * delta1 * img_input';
    W3 = W3 + alpha * delta2 * y1';
end